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On the signature of nearby superclusters and voids in the 
Integrated Sachs- Wolfe effect 
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ABSTRACT 

Through a large ensemble of Gaussian simulations and a suite of large-volume N- 
body simulations, we show that in a standard ACDM scenario, supervoids and su- 
perclusters in the redshift range z £ [0.4,0.7] should leave a small signature on the 
Integrated Sachs Wolfe (ISW) effect of the order ~ 2/iK. We perform aperture pho- 
tometry on WMAP data, centred on such superstructures identified from SDSS LRG 
d ata, and find amplitud es at the level of 8 - 11 fiK - thus confirming the earlier work 
of iGranett etaH (|2008l) . If we focus on apertures of the size ~ 3.6°, then our simula- 
tions indicate that ACDM is discrepant at the level of ~ 4 a. However, if we combine 
all aperture scales considered, ranging from 1°— 20°, then the discrepancy becomes 
~ 2 a. Full-sky ISW maps generated from our A^-body simulations show that this dis- 
crepancy cannot be alleviated by appealing to Rees-Sciama (RS) mechanisms, since 
their impact on the scales probed by our filters is negligible. We perform a series of 
tests on the WMAP data for systematics. We check for foreground contaminants and 
show that the signal does not display the correct dependence on the aperture size ex- 
pected for a residual foreground tracing the density field. The signal also proves robust 
against rotation tests of the CMB maps, and seems to be spatially associated to the 
angular positions of the supervoids and superclusters. We explore whether the signal 
can be explained by the presence of primordial non-Gaussianities of the local type. 
We show that for models with = ±100, whilst there is a change in the pattern 

of temperature anisotropies, all amplitude shifts are well below < 1/iK. If primordial 
non-Gaussianity were to explain the result, then /^L al would need to be several times 
larger than currently permitted by WMAP constraints. 

Key words: cosmology: observations - cosmic microwave background - large-scale 
structure of the Universe - galaxies: clusters: general 



1 INTRODUCTION 

Over the last fifteen years, evidence has been mounting 
from various cosmological probes to support the case for 
an accelerating universe. It can be argued that the first 
compelling evidence for this arose from the study of the 
light curves of distant type la supernovae (|Riess et al.lll998l ; 
iPerlmutter et alJll999h . Currently, the strongest support for 
this picture comes from the combination of observations 
of the Cosmic Microwave Background radiation (hereafter 
CMB), and from measurements of the clustering of galaxies. 
From the CMB side, the Wilkinson Mi crowave Anisotropy 
Probe (hereafter, WMAP) experiment (|Spergel et al.ll2003l . 



120071 : iKomatsu et"all l201lf FI provided a precise measure- 
ment of the angular size of the sound horizon at recombina- 
tion, which supported the case for a spatially flat universe. 
On the clustering side, data from surve ys like the 2-degree 
Field Redshift Survey (|Cole et al.ll2005n and the Sloan Dig- 
ital Sky Survey (hereafter SDSS) (|Tegmark et~ai]|2006F) re- 
quired the density in matter to be sub-critical, hence leading 
to the inference that ~70% of the current energy density of 
the Universe is in a form of energy that behaves like a cos- 
mological constant, and so acts as a repulsive gravitational 
force. The energy density driving the accelerated expansion 
is unknown and so has been dubbed Dark Energy (hereafter 
DE). Uncovering the true physical nature of the DE is one 
of the main targets for many ongoing and upcoming surveys 
of the Universe. 
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Standard linear cosmological theory states that if the 
Universe undergoes a late-time phase of accelerated expan- 
sion, then gravitational potential wells on very large scales 
(> 100 h Mpc) will decay. This evolution of the potential 
wells introduces a gravitational blueshift in the photons 
of the CMB that is known as the Integrated Sachs- Wolfe 
effect (ISW). This effect constitutes an alternative window 
to DE, and can be directly measured by cross-correlating 
CMB maps with a set of tracers for the density fi eld, which 
sources the potentials (|Crittenden fc TurokHl996T ). As soon 
as the first data sets from WMAP were released, several 
works claimed detections of the ISW at various levels 
of significance (Scranton et al. 2003; Fosalba e t alJ 120031 : 
Boughn fc Crittenden! 12004 iFosalba fc Gaztafiagal 12004 



2004 



Nolta et al 
20051; ICabre et al l 



Afshprdi et al 
l200d : 



2004 



Padmanabhan et al.l 
Giannantonio et al.l 120061 ; 



Vielva et all 120061 ; iMcEwen et~ai1 120071 ). Subse- 



quent analysis has led some researchers to be more 
cautious about interpreting t hese earl y detections 
llHernandez-Monteagudo et al. 120061; iRassat et al.l 



20071; iBieibv et all l2010l: lLopez-Corredoira et all l20ld ; 

Hcrnandez-Montcagudo 2010; Francis fc Peacock! |2010D . As 
emphasized by IHernandez- Monteagudo! (|2008l ), the true 
ISW effect should only be detectable for deep galaxy surveys 
that cover a substantial fraction of the sky. However an er- 
roneous interpretation of ISW cross correlation studies may 
be obtained from systematic errors, such as residual point 
source emission in CMB maps or pr esence of spurious power 
on la rge angular scales (see also IHernandez- Monteagudo! 
|2010| . for details on the WMAP - NVSS cross-correlation 
analysis.). In particular, the issue of excess power on 
large scales has been noted in several works llHo et al.l 
20081; iHernandez-Monteagudol l20ld ; iThomas et all l201ll ; 



Giannantonio et al.ll2012l ). but it is not yet fully accounted 



for. 

Amongst the subsequent ISW cross-correlation studies 
in the literatu r e, the re is the particularly puzzling work of 
IGranett et al.l (2008, hereafter G08). Their analysis yielded 
one of the highest detection significances in the literature. 
G08 implemented the following novel approach: they pro- 
duced a catalogue of superclusters and supervoids from 
SDSS data, and stacked WMAP-filtered data on the po- 
sitions of these structures with apertures of the order ~ 4° . 
They obtained a ~ 4 a ISW detection. Unlike previous works 
on the subject, the analysis focused on a particular subset 
of the available large-scale structure (hereafter LSS) data. 
Subsequent studies have investigated the origin of the sig- 
nal and assessed its compatibility with the ACDM scenario, 
dPapai fc Szapudil |2010| ; iPapai et all l201ll ; iNadathur et all 
120121 ). While some works found the G08 results compatible, 
other s found the measured amplitudes too high to be consis- 
tent (|Granett et al.ll2009l ; INadathur et al.ll2012l ). Currently, 
the G08 results remain unexplained. 

In this work we shall attempt to shed new light on this 
problem. The paper breaks down as follows: In [J2]we give a 
brief theoretical overview of the ISW effect, underlining the 
expectation for the ISW signal from a cross-correlation anal- 
ysis of the data used in G08. In iJ3]we perform a cross-check 
of the G08 results. In [JJwe test whether the G08 results are 
consistent with expectations for the ACDM model. This is 
achieved by using a large ensemble of Gaussian Monte-Carlo 
simulations. We also generate full-sky nonlinear ISW maps 



by ray-tracing through a suite of A^-body simulations. In SjS] 
we determine the level of significance at which the results 
disagree with the ACDM paradigm. We explore systematic 
errors in the foreground subtraction for WMAP. We also 
investigate whether the excess signal is consistent with pri- 
mordial non-Gaussinities of the local type. Finally, in [J6]we 
summarize our findings and conclude. 

Unless stated otherwise, we employ a reference cos- 
molo gical model consistent with WMAP7 (Ko matsu et al.l 
2011): the energy-density parameters for baryons, CDM and 
cosmological constant are fib = 0.0456, fi c dm = 0.227, 
Qa = 0.7274; the reduced Hubble rate is h — 0.704; the 
scalar spectral index is ns = 0.963; the rms of relative mat- 
ter fluctuations in spheres of 8 h~ l Mpc radius is as = 0.809, 
and the optical depth to last scattering is r = 0.087. 



2 THEORETICAL PERSPECTIVE 
2.1 The ISW effect 

Observed CMB photons are imprinted with two sets of 
fluctuations: primary anisotropies, sourced by fluctuations 
at the last-scattering surface, and secondary anisotropies, 
induced as the photon propagates through the late-time 
clumpy Universe. The physics of the primary anis otropies 
is well understood (|Dodelsonl 120031 ; IWeinbergj I2008T ) . There 
are a number of physical mechanisms that give rise to 
the generation of secondary an isotropies (for a review see 
iThe Planck CollaboratiorJl2006l ) and one of these is the red- 
shifting of the photons as they pass through evolving gravi- 
tational potentials. The linear vers ion of this effect is te rmed 
the Integrated Sachs- Wolfe effect (|Sachs fc Wolfdl 19671 ) and 
its nonlinear counter- part is termed the Rees-Sciama effect 
(|Rees fc Sciamall 19681 ). 

The observed temperature fluctua tion induced by grav i- 
tational redshifting may be written as (|Sachs fc Wolfelll967l ): 



AT(n) 

To 



(1) 



where n is a unit direction vector on the sphere, $ is the 
dimensionless metric perturbation in the Newtonian gauge, 
which reduces to the usual gravitational potential on small 
scales, the 'over dot' denotes a partial derivative with respect 
to the coordinate time t from the FLRW metric, \ is the 
comoving radial geodesic distance \ = J cdt/a(t), and so 
may equivalently parameterize time. The symbols to and 
iis denote the time at which the photons are received and 
emitted, i.e. the present time and last scattering, c is the 
speed of light and a(t) is the dimensionless scale factor. 

On scales smaller than the horizon, relevant to our simu- 
lation boxes, the perturbed Einstein equations in Newtonian 
gauge lead to a perturbed Poisson equation . This enables u s 
to relate potential and matter fluctuations |Dodclson 20031 ) : 



V 2 $(x;t) = 4irGp(t)S(x;t)a 2 (t) 



(2) 



where p(t) is the mean matter density in the Universe 
and the density fluctuation is defined S(x.;t) = [p(x, i) — 
p(t)]/p(t). This equation may most easily be solved in 
Fourier space: 



*(k;t) = -4nGp(t)a 2 (t)^fi 
k 



(•3) 
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Differentiation of the above expression gives 



6(k;t) = ^tt nM Hik- 



H(t) 
o(t) 



J(k;i) 



<$(M) 



o(t) 



(4) 



where we used the fact that [a 3 (t) p(t)] is a time- 
independent quantity in the matter-dominated epoch 
and r2 m = Sl c dm + fib • In the above, we also defined 
H(t) = a(t)/a(t) and Q m (t) = p(t) / p clit (t), with p CTit (t) = 
3H 2 (t)/8nG. All quantities with a subscript are to be eval- 
uated at the present epoch. In the linear regime, density 
perturbations scale as 5(k, a) = D(a)8(k,ao). Inserting this 
relation into Eq. ([5]) gives, 



*(k ; t) = -a 



noHgk 



>. H(t) 
o(t) 



dlogD 



dloe 



<S(k;a) 



(5) 



In the matter-dominated phase, the Universe expands as in 
the Einstein-de Sitter case and consequently density pertur- 
bations scale as D(a) oc a. Thus, for most of the evolution 
of the late-time Universe the bracket, [1 — dlog D/dlog a] is 
close to zero. In the ACDM model it is only at relatively 
late times that this term is non-zero. 

Alternatively, in the nonlinear regime S(k,a) 7^ 
D(a)8(k,ao), and this gives rise to additional source s for 
the heating and cooling of photons (|Smith et alj |2009). 

2.2 Expectations for cross-correlation analysis 

The expected signal for the ISW-density tracer cross- 
correlation analysis has be e n described by several authors 
JCrittenden fc Turokl 1 19961 ; Irlernandez-Monteagudd |200S| ; 
ISmith et al.ll2009l ) and here we simply quote the main re- 
sults. For a survey of galaxies (or more generically, objects 
that are biased density tracers) the ISW-density angular 
cross-power spectrum can be written: 



C. 



g.ISW 



k dk P m {k) 



■ I d X ijt{kxi)xiW g (xi)n g (xi)Hxi,k)D( X i) 
h 

v Hi , , -3fi m # 2 d(D(x2)/a(x2)) (R) 
x J Q dxiJe(k X 2)— ^ ,(6) 

where P m (k) denotes the linear matter power spectrum at 
the present time, the super-script g refers to the objects 
(galaxies in most cases, but in ours it will refer to voids or su- 
perclusters) probing the gravitational potential wells. In the 
above n g (x) denotes the average comoving object number 
density; W g (x) denotes the instrumental window function 
providing the sensitivity of the instrument to the objects at 
distance x\ the je(x) are the usual spherical Bessel functions 
of order I; and b(x, k) denotes the bias of the tracer popula- 
tion with respect to the matter density field. Th is last factor 
may be both a function of time x an d scale k (|Smith et al.l 
120071 ). but for simplicity we shall simply assume that the 
bias equals unity for the objects in our catalogue. 

In this measurement the primary CMB signal acts as 
noise. If we take this into account and also the variance due 
to the density field, then the S/M with which we expect to 
detect the cross-correlation for a given multipole is: 



(cf ISW ) 2 x (2e + i)f sky 
C? MB (Cf + l/fig) + (cf ISW )' 



where Cf MB and Cf represent the auto power spectra of 
the temperature fluctuations and the density tracers; n g de- 
notes the average angular object number density and hence 
the term l/n g denotes the shot-noise. In our forecast we 
shall neglect this term (hence the predicted significance will 
be slightly over optimistic). The factor f s k y refers to the 
fraction of the sky jointly covered by both the CMB and 
the object surveys. The cumulative S/M for all harmonics 
smaller than I can then be written: 



IT**- 



(8) 



where the monopole and dipole have not been included. 



2.3 Expectations for ISW from data used in 
Granett et al. 

We now turn to the question of what S/M should G08 
have expected to find in the their data. G08 used a 
sample of 1.1 million Luminous Red Galaxies (hereafter 
LRG) selected from the SPSS data release 6 (hereafter 
DR6) l|Adelman-McCarthv et~afl 120081 ). This sample cov- 
ered roughly 7500 square degrees and spanned a redshift 
range of (0.4 < z < 0.75). From this sample they identified 
regions as supervoids or superclusters using the algorithms 
Z0B0V and V0B0Z 0, respectively. The significance of these 
regions was chosen to be at least at the 2cr-level relative to a 
Poisson sample of points. G08 selected the largest 50 super- 
clusters and supervoids for their analysis. Their catalogue is 
publically availably. 

In the left panel of Figure [1] the dot-dashed and dot- 
ted histograms denote the redshift distributions of super- 
voids and superclusters, respectively. The thick solid line 
represents an analytic fit that we have constructed, which 
attempts to be a compromise between the two. Our fitting 
function extends to slightly lower redshifts than the G08 
catalogue, and this will translate into a slightly higher ISW 
prediction. 

The middle panel of Figure Q] presents the angular 
cross-power spectrum as predicted by linear theory for the 
void/supercluster catalogue (solid line). Recall that we are 
taking b — 1 for both voids and clusters. In reality the 
void/cluster regions will be anti-biased/biased and so the 
two signals will differ. However, since here we are more con- 
cerned with the S/N this simplification does not matter, 
especially since we are neglecting the effects of shot-noise 
on the cross-spectra. We also point out that our predictions 
for the ISW-induced cross correlation signal of superclusters 
and supervoids in the SDSS sample are quite similar to the 
predictions for the cross-correlation of AGN s (z < 2) in the 
NVSS catalogue (after adopting the model of lHo etal 1 120081 . 
see the dashed line in the middle panel of Figure [TJ . 

The right panel of Figure [1] presents the predictions for 
the cumulative S/M for the SDSS supercluster and super- 
void analysis. These predictions (solid black line) show that, 
for a ACDM model, one should expect no more than ~ 1.3<r 
significance. In contrast, the prediction for NVSS (dashed 



http : //skysrv.pha. jhu. edu/~neyrinck/voboz/ 



http : //if a.hawaii . edu/cosmowave/supervoids/ 
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Figure 1. Left panel: Redshift distribution of supervoids and superclusters for the catalogue provided by G08 from the SDSS DR6 
data. The dotted histogram denotes supervoids and the dot-dashed histogram denotes superclusters. The solid line represents an analytic 
approximation, which is a compromise between both histograms. Middle panel: angular cross-power spectrum between ISW and the 
projected density. The solid line denotes the predictions for both supervoids and superclusters, modulo the bias being unity. The dashed 
curve represents what one would expect for QSO's in the NVSS survey. Right panel: Cumulative S/Af below a given multipole I. The 
solid line denotes the prediction for matter within the redshift distribution sampled by the G08 catalogue. The dashed line presents the 
predictions for QSO's in NVSS. 



line) is close to ~5.5 (obtained after also neglecting shot 
noise). This is not surprising, since the G08 supervoid and 
supercluster catalogue is relatively shallow, spanning the 
redshift range (0.4 < z < 0.7), and covers only a modest 
fraction of the sky (/ s k y — 0.18). The NVSS is instead sig- 
nificantly deeper and wider. We hence conclude that had 
G08 applied a standard ISW cross-correlation analysis to 
their data, then in the framework of the ACDM model, there 
would have been very little chance for detecting any genuine 
signal at high significance. 



3 CROSS-CHECKING GRANETT ET AL. 
3.1 CMB data 

The WMAP experiment scanned the CMB sky from 2001 
until 2010 in five different frequencies, ranging from 23 GHz 
up to 94 GHz. The angular resolution in each band improves 
with the frequency, but it remains better than one degree 
in all b ands. The S/Af is greater than one for multipoles 
£ < 919 |jarosik et alj|20ld ). and in particular, on the large 
scales of interest for ISW studies, the galactic and extra- 
galactic foreground resid uals are below th e 15pK level out- 
side the masked regions (Gol d et al.ll201lh . 

We concentrate our analysis on the foreground-cleaned 
maps corresponding to bands Q (41GHz), V (61GHz) and W 
(94GHz), after applying the conservative foreground mask 
KQ75y7, which excludes ~25% of the sky. At the scales of 
interest, instrumental noise lies well below cosmic variance 
and foreground residuals, and hence will not be considered 
any further. The ISW is a thermal signal whose signature 
should not depend upon frequency and hence should remain 
constant in the three frequency channels. All of the WMAP 



data employed in this analysis were downloaded from the 
LAMBDA site! 



3.2 Supercluster and Supervoid data 

For our tracers of the LSS, we use the same supercluster and 
supervoid catalogue as used by G08. As described earlier, 
the catalogue was constructed after applying the ZOBOV and 
V0B0Z algorithms to search for supervoids and superclusters 
in the LRG sample extracted from SDSS DR6, respectively. 
G08 used the 50 largest supervoids and superclusters. They 
claimed that this cut yielded the highest statistical signifi- 
cance, in that it minimized the contamination from spurious 
objects, whilst at the same time it provided sufficient sam- 
pling to beat down the intrinsic CMB noise. 



3.3 Methodology 

In their approach G08 have applied a top-hat compensated 
filter or Aperture Photometry (AP) method to the CMB 
map(s) positions of voids and superclusters. This filter sub- 
tracts the average temperature inside a ring from the aver- 
age temperature within the circle limited by the inner radius 
of the ring. In order to have equal areas in both cases, the 
choice of the outer radius of the ring is y2R, with R the 
inner radius of the ring. In this way, fluctuations of typical 
size R are enhanced against fluctuations at scales smaller or 
larger than such radius. Although G08 present results for 
apertures ranging from 3° up to 5°, most of the conclusions 
are driven from the R = 4° choice, for which highest sta- 
tistical significance is achieved: they find that AP stacks on 
the position of voids (superclusters) yield a decrement (in- 
crement) of ~ — 11.3/iK (7.9 j«K) at 3.7 (2.6) a significance 



4 http://lambda.gsfc.nasa.gov 
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Figure 2. AP filter outputs versus aperture size. Red, green and 
blue symbols refer to WMAP's Q, V and W bands, respectively. 
Squares and circles refer to superclusters and voids, respectively. 

level. However, it turns out that, according to G08, the typ- 
ical size of clusters and voids are ~ 0.5° and 2°, respectively, 
which seem to lie at odds with the aperture choice of 4° . Po- 
tentials are known to extend to larger scales than densities, 
and it is a priori unclear which aperture radius should be 
used. This fact motivates a systematic study in a relatively 
wide range of aperture radii. 



3.4 Stacking analysis on real data 

We apply the G08 method on the Q, V, and W bands of 
the WMAP data using the SDSS supercluster and super- 
void catalogues. We have considered AP filters in 15 log- 
arithmically spaced bins in the angular range 1°— 20°. The 
filters were placed on the centers of the objects, as they are 
provided by the catalogue. 

Figure [2] displays the results for the stacked signal as a 
function of the AP filter aperture size in degrees. The red, 
green and blue symbols refer to results from the Q, V and W 
bands, respectively. We clearly see that there is practically 
no frequency dependence. The error bars are computed af- 
ter repeating the analysis on 30 random sets of 50 objects 
placed in the un-masked region of the sky. Our findings are 
in good agreement with those of G08: voids and superclus- 
ter regions yield a slightly asymmetric pattern, with voids 
rendering amplitudes of ~ — 11/iK for apertures of ~ 3.6°, 
and superclusters giving rise to increments of ~ 9/iK at 
that same scale. In these two cases, the significance is about 
-3.3 a and 2.3 a, which on combination yields a combined 
significance ~ 4a. Again, this is in good agreement with 
G08. Furthermore, the scale showing highest S/AT is < 4°, 
as the significance rapidly drops for smaller and larger aper- 
tures. Intuitively, this seems to be in contradiction with the 
idea of ISW fluctuations being large-scale anisotropies, since 



Figure 3. AP filter outputs for the Gaussian simulations of the V 
band of WMAP. Results for simulations of other bands are virtu- 
ally identical to these. Black, green and blue symbols and curves 
refer to density threshold choices of v = 2.5, 2.8 and 3.2, respec- 
tively. Circles correspond to voids and squares to superclusters, 
while solid lines display the S/AT at each scale. 

in such case one would expect to attain high S/Af also for 
moderately large (> 5°— 10°) apertures. 



4 ARE THE RESULTS CONSISTENT WITH 
THE VANILLA ACDM UNIVERSE? 

4.1 Gaussian simulations 

What is not clear from the analysis of the previous section, 
is whether the ~ 4a detection from the stacked superclus- 
ter/supervoid regions is consiste nt with what one ex pects 
from the standard ACDM model (|Komatsu et al.ll201lh . We 
now attempt to understand the dependence of the expected 
signal and its errors on the filter size. To that end, we re- 
peat the above analysis on a set of Gaussian simulations of 
both the LSS distribution and the corresponding CMB tem- 
perature anisotropy distribution which would result in our 
ACDM model. 

The CMB maps are constructed in a two-step process: 
first, we generate a Gaussian map of projected density, fol- 
lowing the angular power spectrum built upon the redshift 
window function W s (r) displayed in the left panel of Fig- 
ure [T] This density map is used for (i) constructing a super- 
cluster and supervoid catalogue (see below), and (ii) gen- 
erating an ISW component. This ISW component is cor- 
rela ted to the density map as predicted by Eq. El (see, 
e.g. JCabre et al . 2006; H ernandez-Monteagudd[2008l .for de- 
tails.) Second, we generate the primary anisotropy signal 
at z ~ 1 050. This is taken to be completely uncorrelated 
with respect to the projected density map. This is not ex- 
actly true due to the lensing of the CMB, but it is still a 
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Integrated Sachs-Wolfe effect between z = 0. 0—1.0 




-34.8 



Over density map between z = 0. 366-0. 270: Gaussian smoothing:= 16.9 arcmins 




Figure 4. Top panel: Full-sky, nonlinear ISW map arising from structures between z = 0.0 and z = 1.0, estimated from one of the 
zHORIZON simulations. The sky map was pixelated using the HEALPix package with a resolution of N s ^ c = 256, which corresponds to 
786,432 pixels. Bottom panel: full-sky projected over-density map for the distribution of cold dark matter structures located in the 
redshift shell z = [0.260,0.366] for the same ./V-body simulation. The map was pixelated at the same resolution as the ISW map and 
smoothed with a Gaussian filter of scale FWHM = 16.9 arcmins (roughly the angular size subtended by lh _1 Mpc at z = 0.3). In both 
panels the graticule scale indicates 30° divisions. 



very good approximation on the large angular scales of in- 
terest in this study. The simulated ISW map is then directly 
co-added to the primary CMB temperature map. We have 
checked that the cross-correlations of the simulated density 
and CMB maps are in direct agreement with a numerical 
evaluation of Eq. [6] 

All of the simulated sky maps were generated using the 
equal area pixelization strategy provided by HEALPijJj. We 
take the pixel scale to be ~ 15', corresponding to a HEALPix 
resolution parameter of iV s idc = 256. For each simulated 



LSS map, we smooth the map with a Gaussian aperture 
of FWHM~2° , and identify those peaks and troughs which 
exceed a given threshold v\a\ as being associated with su- 
pervoids and superclusters, where a in this context refers 
to the density field rms. The threshold v is chosen to have 
an object density similar to that of superclusters and voids 
in the real catalogue under the SDSS footprint. To check 
the dependence on the threshold choice, we bracket the pre- 
ferred value of v with two other values, one above and one 
below it. The final choice for the v value set was 2.5, 2.8 and 
3.2. 



http : //www.healpix . org 



To each simulated CMB sky, we also add a noise real- 
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ization following the anisotropic noise model provided in the 
LAMBDA site. We then exclude pixels in accordance with the 
intersection of the WMAP KQ75y7 sky mask and the SDSS 
DR6 data footprint. 

Figure [3] presents the ensemble-averaged results ob- 
tained from 1000 Gaussian Monte-Carlo simulations. As ex- 
pected, for the three adopted thresholds, overdensities (or 
positive excursions in the projected density map) yield a 
positive signal for the stacked aperture analysis (as displayed 
by the squares in the plots), whereas underdensities yield 
negative ones (circles in the plot). The Monte-Carlo simula- 
tions also enable us to compute the variance on the measure- 
ments. On taking the ratio of the mean signal and the rms 
noise, we obtain direct estimates of the S/M for each given 
aperture bin. The coloured solid lines in the plot present 
our direct measurements of the S/M. Thus we clearly see 
that the scatter induced by the CMB generated at z ~ 1050 
is the dominant source of noise, keeping the S/M for each 
angular bin below unity. 

For lower thresholds there is more area covered and in- 
tuitively one would expect a higher S/M, as it seems to 
be the case. Our simulations also provide higher ISW am- 
plitudes for higher thresholds, and this makes sense since 
deeper voids/potential wells should have a stronger impact 
on CMB photons. Nevertheless, in all cases typical ampli- 
tudes remain at the level of 1-2 /iK. The aperture at which 
the AP outputs provide the highest amplitude does not show 
any strong dependence upon the threshold u, and seems to 
lie in a wide angle range within [3°, 8°]. For apertures larger 
than 10°, the S/M for the lowest threshold starts dropping 
slowly, and becomes half of its maximum value at an aper- 
ture of 20°. For higher thresholds this decrease is found to 
be even shallower. 

4.2 Generation of nonlinear ISW and density 
maps from A-body simulations 

The previous subsection has shown that the Gaussian re- 
alizations of the ACDM universe are in tension with the 
excess signal found by G08. One weak point in the above 
analysis is that the density and late-time potential field are 
not necessarily well described by a Gaussian process, since 
nonlinear evolution under gravity drives the initially Gaus- 
sian distribution of density fluctuations towards one that is 
non-Gaussian at late times. In order to test whether nonlin- 
ear evolution could explain the excess signal seen by G08, 
we now turn to the challenge of constructing fully nonlin- 
ear maps of the density field and ISW effect from TV-body 
simulations. 

The 8 simulations that we employ for this task are 
a sub-set of the zHORIZO M simu lations. These simulations 
were used in ISmith et al l £2009) to calculate the expected 
ISW-cluster cross-power spectra. In brief, each simulation 
follows the gravitational evolution of N = 750 3 dark mat- 
ter particles in a box of comoving size L = 1500/i _1 Mpc. 
The cosmological model employed was a flat ACDM model: 
Qmo = 0.25; cr 8 = 0.8; n s = 1.0; h = 0.72, fi b ,o = 0.04. The 
transfer functi on for the simulations was ge nerated using the 
cmbfast code (Scliak & Zaldarriaga 199y). The initial con- 
ditions \^ere_\mn^mm^tied^^h_z_= 49 u sing the code 
2LPT (|Scoccimarrd 1 19981 ; ICrocce et all 120061 ). Each initial 
condition was integrated forward using the publicly available 



cosmological A^-body code Gadget-2 |Springelll2005l ). Snap- 
shots of the phase space were captured at 11 logarithmically- 
spaced intervals between a = 0.5 and a — 1.0. 

In order to generate full-sky nonlin ear ISW maps w e 
roughly follow the strategy described in ICai et all |2010h , 
but with some minor changes. Full details of how we con- 
struct our maps can be found in Appendix fS] In summary, 
we used the density and divergence of momentum fields to 
solve for $ for each snapshot. We then constructed a back- 
ward light-cone from z — 0.0 to z = 1.0 for 4>. We then 
pixelated the sphere using the HEALPix equal-area decompo- 
sition, taking the pixel resolution to be A s id c = 256, which 
corresponds to 786,432 pixels on the sphere. For each pixel 
location, we then fired a ray through the past light-cone of $ 
and accumulated the line-of-sight integral given by Eq. fT]). 
Note that we only consider the ISW signal coming from 
2 < 1, since we do not expect a significant cross-correlation 
between the relatively low-redshift density slices for SDSS 
and the ISW from z > 1. The top panel of Figure|4]shows one 
of the ISW maps that we have generated from the zHORIZON 
simulations. 

We next generated the projected density maps. These 
were done by first constructing the projected density map 
associated with each snapshot at. The density field for a 
given snapshot was obtained as follows. To each snapshot a; 
we associate a specific comoving shell [Xi-i/2, Xi+1/2] (see 
Appendix [A] for more details). We then select all of the 
particles that fall into the shell for that epoch, i.e. the ith 
dark matter particle in the ajth snapshot, is accepted in the 
shell if Xi-i/2 < |xi - Xo| < Xl+1/2, where xj and x are 
the coordinates of the ith dark matter particle and the ob- 
server, respectively. Note that if a given value of x is larger 
than L/2, 3L/2, 5L/2, . . . , then we apply periodic bound- 
ary conditions to produce replications of the cube to larger 
distances. If the particle is accepted, then we compute the 
angular coordinates (8, <j>) for the particle, relative to the 
observer. Given these angular coordinates, we then find the 
associated HEALPix pixel and increment the counts in that 
pixel. The bottom panel of Figure [4] shows the projected 
overdensity map for a thin redshift shell centred on z — 0.3. 
We note that it is hard, by eye, to note any apparent corre- 
spondence between the overdensity and temperature maps. 

At the end we have 11 density maps between 2 = 1.0 
and z = 0.0 that form concentric shells around the observer. 
These shells were then co-added using the weights given by 
our analytic fit to the redshift distribution of superclusters 
and supervoids (recall the solid line in the left panel of Fig- 
ure!]]). The resulting co-added all-sky density maps are then 
smoothed from A^ide = 256 to N s id c = 32, and the positions 
of the 2n, n and n/2 most under- and over-dense pixels 
on this map are recorded. The number n corresponds to a 
number density of extrema that is identical to that of real 
voids and superclusters under the footprint of SDSS DR6. 
Each of those extreme pixels on the A s id c = 32 map is then 
projected back to the A^dc = 256 map on a subset of 64 
higher-resolution pixels, out of which the position of the 
most under- or over-dense pixel is used as the target of the 
AP filter. 
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Figure 5. (Top) Comparison of the average angular power spec- 
tra obtained from the zHORIZON-derived ISW maps (black cir- 
cles) and the theoretical expectation (red solid line). The lack of 
low-A; modes causes a low bias in the low multipole range, and 
non-linear evolution introduces some visible power excess on the 
small scales. (Bottom) Comparison of the average angular power 
spectra from the zHORIZON-derived density contrast maps (black 
circles) with the linear prediction (red solid line), for the redshift 
range z € [0.1, 1]. The agreement with theory is again good on 
intermediate and small scales, but there seems to be again some 
slight low bias on the largest scales/lowest multipoles. 



4.3 Validation tests of zHORIZON derived maps 

Before we apply the analysis methods of G08 to our ISW and 
density maps, we first test the consistency of the maps them- 
selves. To do this, we compute the angular auto-power spec- 
tra of the ISW temperature maps for each of the 8 zHORIZON 
runs. The top panel of Figure [S] presents the ensemble- 
averaged temperature power spectrum for the ISW effect. 
The measurements from the simulations are represented by 
the solid black points. The prediction from linear theory 
is given by the solid red line. For multipoles in the range 
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Figure 6. AP filter outputs for ISW maps derived from the 
zHORIZON simulations. Black, green and blue symbols refer to 
density excursions in the projected density maps having twice, 
same and half the angular number density as voids and superclus- 
ters from G08. Circles correspond to under-densities and squares 
to over-densities. These excursions were identified in projected 
density maps smoothed down to pixels of ~ 2 deg on a side 
(N sido = 32). 



5 < £ < 70, the agreement between the two is excellent. At 
low multipoles {£ < 5) the absence of power in the simu- 
lation on scales larger than L induces a low bias. Instead, 
on scales £ > 70, the non-linear evolution of potentials sub- 
stantially boosts the signal relative to linear by means of the 
Rees-Sciama effect. 

We next compute the angular auto-power spectrum of 
the projected density contrast maps, with the projection ex- 
tending from z = 0.1 up to z — 1 (bottom panel of Figure[SJ). 
In this case, the projected density angular power spectrum 
shows good agreement with the linear-theory expectations 
in the intermediate- and high- multipole range (5 < £ < 200), 
and some hints for power deficit on the large scales/low £-s, 
which would probably be due to the lack of k modes beyond 
the box size of the simulations. 



4.4 Aperture analysis of the zHORIZON maps 

Having validated the simulated maps, we now repeat the 
AP analysis. Since here we have both full-sky density and 
temperature maps, we prefer not to apply any sky mask. 
Thus, these predictions will not be affected by incomplete 
sky coverage. We repeat the steps described in fj3]for finding 
the locations of the density peaks/troughs in each of the 
simulated smoothed maps. Then, as before, we apply the 
AP filters to the selected centroid positions for the 8 ISW 
maps from the zHORIZON simulations. 

Figure HJ] presents the results from this analysis. The 
square and circular symbols denote the results for our ef- 
fective supercluster and supervoid regions, respectively. The 
blue, green and black colours correspond to the set of ex- 
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treme pixels which have half, equivalent and double the an- 
gular number density of the real supervoids and superclus- 
ters found in G08's analysis. We see that the peak of the 
average AP output has a temperature of the order 2 fiK, 
and occurs for apertures of scale ~ 6-7° . 

On comparison with the predictions from our Gaussian 
simulations, we find that the fully nonlinear ISW maps are 
in close agreement (c.f. H4.ip . The are however small differ- 
ences. The peak signal is shifted to slightly larger scales for 
the full nonlinear case. Also, the shape of the curves obtained 
from the zHORIZON simulations appears smoother than in the 
Gaussian simulation case, which shows dips and troughs that 
are absent in Figure [6] We believe that these small differ- 
ences are likely a consequence of the fact that the Gaussian 
simulations include intrinsic CMB noise and possess a sky- 
mask. However, the most important point to note, is that 
on angular scales of the order 3-4°, the Gaussian and fully 
nonlinear simulations are in close agreement: the difference 
induced by adopting a slightly different cosmological model 
should introduce changes in the ISW amplitude at the 2 % 
level, and the ISW generated beyond z = 1 seems to have 
little impact as well. We thus conclude this section by not- 
ing that that the excess temperature signal found by G08, 
and now confirmed by us in |3] appears to be incompatible 
with the evolution of gravitational potentials in the standard 
ACDM model. 



5 SIGNIFICANCE, SYSTEM ATICS AND 
ALTERNATIVE MODELS 

5.1 Estimating the significance 

The direct comparison of Figures[2]and[3]reveals clear differ- 
ences between the observed data and theoretical predictions. 
Not only is the amplitude of the maximum signal in the real 
data a factor of ~5 times larger than the average in the 
Gaussian simulations, but the dependence of the signal on 
the filter scale shows a different shape. More quantitatively, 
the results from the WMAP data for an aperture size of 
3.6°, are of the order ~ 3.3 a away from the supervoid sim- 
ulation average, and ~ 2.1 a away from the average for the 
case of superclusters. In terms of probability, for an aperture 
scale of 3.6°, only one out of the 1000 simulations possessed 
an ISW signal, from supervoid regions, with a temperature 
decrement lower than the one found in real data, and none 
of the simulations for superclusters exceeded the value ob- 
tained for the real data. Taken at face value, this analysis 
seems to exclude the Gaussian ACDM hypothesis at > 3a 
significance. However, this is an a posteriori estimate, since 
we have neglected the fact that we also looked for a signal 
at other aperture scales. 

If we include the measurements from the 15 different 
angular aperture scales between 1-20° and also if we take 
into account their covariance, then the significance drops. 
Under the assumption of Gaussian statistics, we find that 
the WMAP outputs for voids produce a Xvoids — 23.0 
(?Vof = 15). The corresponding figure for superclusters is 
Xsuperciustcrs — 24.0 (n do f = 15). If we treat these two con- 
straints as being independent then their combination yields 
Xboth = 46.9 (ndof = 30). In terms of probability, this means 
that the WMAP data has a 0.025 probability (i.e. < 3% 
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Figure 7. Rotation test on AP filter outputs for an aperture size 
of 3.6°. Blue circles (red squares) correspond to voids (superclus- 
ters). 

chance) of being consistent with the evolution of gravita- 
tional potentials in the ACDM model. 

5.2 Tests for systematics 

We next test the possibility of systematics in WMAP data 
giving rise to the observed signal. 

• Rotation test: We first conduct a rotation test as- 
sessing the statistical significance of the AP outputs for a 
3.6°scale. In this way, we probe the possibility of any other 
signal (apart from CMB) contributing to the uncertainty of 
the AP filter outputs and hence modifying the statistical 
significance found for WMAP data. We rotate in galactic 
longitude (in steps of 9°) the AP filter targets with respect 
to the real positions of supervoids and superclusters. In the 
absence of systematics, this should provide AP outputs com- 
patible with zero. 

Figure [7] presents the results from this analysis. The blue 
circular symbols represent the supervoid regions and the 
red square symbols denote the superclusters. At zero ro- 
tation lag we clearly obtain a signal of higher amplitude 
than in any other rotation bin. We have verified that this 
signal does not arise as a consequence of a small subset of 
the supervoids/superclusters. Instead, the signal is approx- 
imately evenly distributed among all structures. From the 
sample of rotated bins only, the estimated significance for 
the 3.6°aperture is 4.1a for voids, 2.7 a for superclusters 
and 3.8 a combined. This is within 1-a from the significance 
levels obtained with the Gaussian simulations. 

• Density-dependent contaminant: Another possible 
systematic is some combination of contaminants present in 
the superclusters that increases the emission in these struc- 
tures in a frequency independent manner. This could happen 
if the CMB cleaning algorithms remove only the frequency 
varying part of the contaminant signal, but leave behind a 
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ISW effect between z = 0. 0—1.0: Difference between FNL= + 100 and Gaussian case 
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Figure 9. Relative difference between full-sky ISW temperature maps for models with primordial non-Gaussianity and Gaussian models 
for flat ACDM universes. Top panel: presents AT ISW (/Jf£ al = +100) - AT ISW (/J,°£ al = 0). Bottom panel: presents the same but for 
AT ISW (/Jj£ al _ „100)-Ar isw (/J5£ al = 0). All maps were smoothed with a Gaussian filter of FWHM = 1.0°, which roughly corresponds 
to the scales at which the signal for the AP filters peaks. The large-scale temperature fluctuations appear slightly cooler/hotter in Universes 
with /jJtJ positive/negative. 



DC level. In such a scenario, one might assume that the sig- 
nal is linearly responding to fluctuations in the projected 
matter density. We incorporate this effect into our Gaussian 
simulations by substituting the total CMB map by a signal 
that is proportional to the projected matter density field 
in our Gaussian simulations. This should provide us with 
some idea of the scale-dependence for a contaminant of this 
nature. 

Figure[8l presents the results from this exercise. The sym- 
bol styles and point colours are the same as in Figure 
The recovered shape resembles that for a standard matter 
auto-correlation function, and is remarkably different from 
that of Figure [2] Actually, the profile from real data in Fig- 



ure [5] seems to be an intermediate case between the sce- 
nario depicted in Figure [8] and the theoretical predictions 
of Figure [3] However, if most of the observed amplitude at 
3.6° (~ 8 — 2 = 6/iK out of the total ~ 8/iK observed) is 
to be caused by this type of contaminant, then the scale de- 
pendence of the output should accordingly be much closer 
to the one shown in Figure [5] and this is not the case. Note 
that we have assigned ~ 2 to ISW in this estimation. If 
contaminants in the position of the AP targets were Poisson 
distributed, then the profiles obtained in Figure [S] would ap- 
proach zero faster as aperture radii increase, yet in stronger 
disagreement with Figure [2] 

• Superstructure selection effects: It might be ar- 
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Figure 8. Pattern induced versus AP filter scale by a contami- 
nant following the local matter density. 



gued that our identification of under- and over-dense re- 
gions in our simulated density maps does not match well 
the selection of supervoids and superclusters in the G08 
data. Whilst in principle both approaches should generate 
the most under- and over-dense regions in volume limited 
samples of the universe, differences in the exact details of 
the Z0B0V and V0B0Z implementations may introduce subtle 
but important differences in the selection of regions. This 
would have an impact in the final definition and selection 
of targets for the AP filters. Properly addressing this issue 
would require a thorough implementation of the Z0B0V and 
V0B0Z algorithms to semi-analytic galaxy samples embedded 
in our simulations (and this goes beyond the scope of the 
current work). Having said this, the stability of our results 
with respect to the actual density peak threshold adopted, 
suggests that this would not critically affect our conclusions. 



5.3 Primordial non-Gaussianity 

In this section we explore whether the G08 result is com- 
patible with a non-Gaussian distribution for the primordial 
potential perturbations. In particular, we consider the well 
known local model for primordial non-Gaussianity, charac- 
terized by a quadratic correction to the Gauge invariant 
Bardeen's potential perturbation fsee lKomatsu et alj|20lj , 
and references therein): 

$ NG (x) = $ G (x) + /^ L cal [$G(x) 2 -($|(x))] , (9) 

where &g(x) is the Gaussian potential perturbation after 
matter radiation equality, scaled in terms of units of c 2 to 
yield a dimensionless quantity. Following standard conven- 
tion, $ NG (x) = -$ Ncwton (x) (i.e. the $ in Eq. ©). The 
term (&Gr{x)) is subtracted to ensure that <E>ng is a mean 
zero field. In linear theory the typical fluctuations are of the 
order $ng ~ 1CP 5 , and so the non-Gaussian corrections are 
of the order ~ 0.1%(/^/100)($G/10 _5 ) a . 

In order to explore the observable consequences of such 



a modification, we have generated a set of simulated ISW 
maps with Gaussian initial conditions, i.e. /^ al = 0, and 
with non-Gaussian initial conditions /Jj°L al = { + 100, —100}. 
These maps were generated from A-body simulations seeded 
with Gaussian and non-Gaussian initial conditions follow- 
ing the methodology of AppendbJXl The simulations tha t 
we employ were fully described in lDesiacques et al.l ((2009). 
In brief, these were performed using Gadget-2, and fol- 
lowed N — 1024 3 dark matter particles in a box of size 
L = 1600 h~ 1 Mpc. The cosmological mode l of the simula- 
tions was consistent with the WMAP5 data (|Komatsu et al.l 
2009). We use a sub-set of these simulations that were used 
elsewhere for gravitation al lensing analysis (|Marian et al.l 
l201ll ; iHilbert et al.ll20rl ). The simulations were set up to 
have the same initial random phases for all three models, 
this enables us to cancel some of the cosmic variance and so 
permit us to better explore the model differences. 

The top panel of Figure [JJ] presents the differences 

between the ISW temperature maps in a Universe with 
/IjT' _ +100 and jiocai _ q The bottom panel shows 

the same but for the case = -100 and /^ al = 0. 

Note that all of the maps were smoothed with a Gaussian 
filter of FWHM = 1° before being differenced. Note also 
that we have only included the ISW contributions between 
z = 0.0 and z = 1.0. Clearly, the presence of primordial 
non-Gaussianities can induce shifts in the temperatures of 
the peaks and troughs of the distribution. However, these 
shifts are modest for /Jj 3 ' = ±100, leading to changes that 
are < l^iK. This suggests that values of /J(™ al on the order 
of ~ 1000, might be able to explain the AP analysis of the 
WMAP results. However, such large values for /Jj ™ 1 would 
be grossly inconsistent with the values of /nl obtained 
from the CMB temperature bispec trum, which currentl y 
gives -10 < /N°L al < 74 (95% C.L., iKomatsu etafl l201ll ). 
It therefore seems unlikely that the scale-independent local 
model of primordial non-Gaussianity is the correct explana- 
tion for the excess signal. 



6 CONCLUSIONS 

In this work we have studied the imprint of superclusters 
and supervoids in the temperature map of the CMB from the 
WMAP experi ment. Our work f urther explores the signature 
first detected in lGranett et all (|200Sl . G08). 

In Sj3]we theoretically showed that if G08 had applied 
a standard angular cross-power spectrum analysis of the su- 
perstructures they found in the SDSS LRG data, then the 
expected significance for a ACDM model should have been 
< 1.5 a 

In |3] we cross-checked the G08 analysis directly and 
found identical conclusions: on scales ~ 3.6° there was a 
~ 4<7 detection significance for excess signal associated with 
the supervoids and superclusters. 

In 33] we performed a series of tests exploring whether 
these findings are consistent with the standard ACDM 
model. Gaussian Monte-Carlo simulations of the ISW ef- 
fect and the LSS were unable to produce such large signals. 
We then investigated whether this was a consequence of our 
simplified Gaussian simulations. We did this by generating 
fully non-linear ISW maps from large volume Af-body simu- 
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lations. These simulated maps confirmed the findings of the 
simpler Gaussian simulations. 

In Sj5]we used the Gaussian Monte-Carlo simulations to 
explore the significance of the deviations from the ACDM 
model found in the WMAP data. We found that for aperture 
photometry analysis of the maps on scales 3.6°, only 1 out of 
the 1000 mocks had a supervoid signature compatible with 
the G08 result, while none of the supercluster regions did. 
On taking into account the 15 aperture scales examined, the 
significance of the discrepancy dropped to < 3% chance of 
the result being consistent with the ACDM model. 

We investigated whether the effect could be caused by 
a systematic error in the cleaning of foregrounds in the 
WMAP data. We found that if the signal were to be caused 
by an approximately frequency-independent emission trac- 
ing the density field, then the resulting angular dependence 
would be very different to the measured shape found in the 
WMAP data. 

We explored whether the signal could be generated 
by primordial non-Gaussianities. We considered the local 
model, characterized by a quadratic correction to the pri- 
mordial potential perturbations, with the coupling parame- 
ter /nT'- We found that, for positive/negative, asym- 
metric shifts in ISW temperature maps arise. However, for 
the values of /nl 3 ' = ±100, the changes were < 1/iK (after 
smoothing the maps down to degree scales). Thus values of 
/nl 3 ' an order of magnitude higher would be required to ex- 
plain the G08 result, and they would be clearly inconsistent 
with current constraints on /NL al from WMAP. 

It is possible that the G08 result may also be explained 
by other more exotic scenarios, e.g., non-Gaussianity arising 
from the presence of a non-zero primordial equilateral or or- 
thogonal model bispectrum (a consequence of non-standard 
inflationary mechanisms); alternatively it might arise as a di- 
rect consequence of modifications to E instein's general the- 
ory of relativity, l|Jain fe Khou rv 2010). However, more con- 
servative scenarios involving some combination of artifacts 
and/or systematics cannot yet be fully discarded. 

In the future, we will look with interest to the results 
from the Planck satellite as to whether this signal represents 
a data artifact, or whether it constitutes a genuine challenge 
to the ACDM model and a window to new cosmological 
physics. 
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APPENDIX A: FULL SKY ISW MAPS FROM 
iV-BODY SIMULATIONS 

In this section we aim to construct full-sky ISW maps us- 
ing a suite of A r -body si mulations . Our approach is sim- 
ilar to that described in ICai et all (|2010h . but with some 
modifications. To be more precise, we aim to compute the 
line-of-sight integral Eq. |l} , but taking into account the full 
nonlinear evolution of "I>. The steps we take to achieve this 
are described below. 



Al Determining <E> 

In order to obtain $ directly from the iV-body simulations, 
we make use of Eq. ((5J , which tells us that our desired 



quantity can be determined from knowledge about S(k,a) 
and S(k, a). In simulations, measuring S(k, a) is relatively 
straightforward, whereas its time deriv ative is more compli- 
cated. As was shown by Seliak (1996), this latter quantity 
may be obtai ned from the perturbed continuity equation 
l|Peeblesll 19801 ): 



V ■ [1 + <S(x; t)] v p (x; t) = -a(i)5(x; t) 



(Al) 



where v p (x;i) is the proper peculiar velocity field. If we 
define the pseudo-peculiar momentum field to be, 



p(x;t) = [l + <S(x;t)]v p (x;i) , 



(A2) 



then in Fourier space we may solve the continuity equation 
directly to find 



<5(k;t) = ik-p(k;t)/a(t) 
Hence, our final expression becomes, 



$(k;t) = F{k) 



Hit) 
a(t) 



*(k;t)- 



zk ■ p(k; t) 



where we defined the function 



(A3) 



(A4) 



(A5) 



Thus in order to estimate <!>, we simply require estimates of 
both the density field and pseudo-peculiar momentum field 
in Fourier space. 

The dark matter density field in an iV-body simulation 
can be written as a sum of Dirac delta functions, 



p(x) = ^2 miS D (x - X;) 



(A6) 



where N is the number of particles and mi is the mass of 
the Zth particle, and we take all particles to have equal mass. 
The density field averaged on a cubical lattice can then be 
obtained through the convolution, 



1 



d xp(x)W(xij fe — x) 



(A7) 



where xy*, labels the lattice point, W represents the 
dimensionless window function of the mass assignment 
scheme. This window function is normalized such that 
Vw = J d 3 x.'W(x - x'). The filter function W that we adopt 
throughout is the 'cloud-in-c ell' charge assignment scheme 
jHocknev fc Eastwood! 1 19881 ). Hence, our estimate for the 
density fluctuation is 



1 + 5(x ljk 



N 



(A8) 



where Aceii = V^/Vw is the total number of grid cells. 

The pseudo-momentum field may be estimated in a sim- 
ilar fashion. For convenience we write, 



p = [1 + <5(x)] u(x)a(f) 



(A9) 
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where u = v p /a is the comoving peculiar velocity field. The 
particle momentum field is then written as 



t(l + *)u] (x) = ^^5 D (x-x ; ) Ui 



(A10) 



This may be convolved with the mass assignment scheme to 
obtain the mesh averaged quantity 



[(l + *)u] (x ijfc ) 



iV w 



Xijfc - X;) 



(All) 



Thus our estimate for the pseudo-momentum field is 

N N 

p(x,jfc) = a(t)— j^- uiW(x ijk - Xj) . (A12) 

! 

The density Fourier modes we re then estimated usin g 
the publicly available FFTW routines | Johnson fc Frigof 2008), 
and each resulting mode was corrected for the convolution 
with the mass- assignment window function. For the CIC al- 
gorithm this corresponds to the following operation: 



<5 d (k) =tf,(k)/Wbra(k) , 



where 



Wbio(k)= [] 



sin [7rfci/2fcN y ] 



[nki/2k-N y ] 



(A13) 



(A14) 



and where sub-script d and g denote discrete and grid quan- 
tities, and where &N y = irNg/L is t he Nyquist frequency, 
and N s is the number of grid cells l|Hocknev fc Eastwood! 
1 19881 ). 

To obtain the real space <I>(x,i), we solved for <i>(k, i) 
in Fourier space using Eq. (IA4[) . set the unobservable k = 
mode to zero, and inverse transformed back to real space. 



A2 Reconstructing the light-cone for $ 

We now wish to construct the past light-cone for the evolu- 
tion of <&, however we only have a finite number of snapshots 
of the particle phase space from which to reconstruct this. 
It is usually a good idea to space snapshots logarithmically 
in expansion factor, and for simplicity we shall now assume 
that to be true. The light-cone can then be constructed as 
follows: 

• Find 4>(xjj&, ai) for every output I on a periodic cubical 
lattice, using the techniques described in Appendix I All 

• Place an observer at the exact centre of the simulation 
cube, xo, and compute the comoving distances from the 
observer to the expansion factors, a;_i/ 2 , at, and a t+1 ^ 2 , 
and label these distances Xi-1/2, Xi an d Xi+i/2- Here, 



log a, 



±1/2 



log ai ± A log a/2 



with A log a being the logarithmic spacing between two dif- 
ferent expansion factors. The comoving distance from the 
observer at do to expansion factor a is given by: 



X(a) 



cda 
a 2 H{a) 



Hence, the intervals [xi-1/2, Xl+1/2] form a series of concen- 
tric shells centred on the observer. 



• Construct a new lattice for the <£> values on the light 
cone. This is done by associating to each snapshot a; , a spe- 
cific comoving shell [xi-1/2, Xi+1/2], and taking only those 
values for $(xjjfc, d;) that lie within the shell: i.e. if 

Xi-i/2 < l x ijfc - xo| < Xi+1/2 , 

then 4>(xijfc, ai) is accepted onto the new grid. Note that if a 
given value of x is larger than L/2, 3L/2, 5L/2, . . . , then we 
use the periodic boundary conditions to produce replications 
of the cube. 



A3 Computing the ISW line-of-sight integral 

Having constructed the backward lightcone for $, we may 
now compute the line-of-sight integral for the ISW effect 
through Eq. ([T]). In fact we use a slightly different form of 
this equation by transforming the integration variable from 
t to log a: i.e. 



AT(n) 

To 



lnao ,, *(n,x) 
a log a- 



H{a) 



(A15) 



The above expression can then be discretized to give, 

, M 



AT(h) 
To 



-^Aloga, 



1=0 



$(n,aj) 
Hifli) 



(A16) 



If we take A log a to be constant, then this becomes, 
AT (n) 2Alog a A $(n,a<) 



^ H 



{(H) 



(A17) 



In evaluating the above expression, we take the number of in- 
tervals to be as large as desired, b ut always evaluating H(aj) 
exact ly through the expression (|Dodelsonl 120031 ; IWeinbera 

Hooi), 

H 2 (a) = Hq [Q,ao + fimoa -3 — (fi m0 + ^ao — l)a~ 2 ] , 

(A18) 

which is valid in the matter dominated epoch for the ACDM 
model. Whereas, for $(0^), we employ the light cone derived 
form from the previous subsection, but interpolate the value 
off the 3D mesh using a CIC like scheme. That is: 

$(x,a) = (1 - hm)(l - hy)(l - h x )9(Xij, h ) 

+ (1 - hx)(l - h y )hA{-x.i,j,k+i) + ■■■ 

)) , (A19) 

where {h x ,h y ,h z } are the x-, y-, and ^-coordinate separa- 
tions of the evaluation point x and the position vector for 
lattice point k), in units of the inter-lattice separation. 

In making this separation of the evolution of H(a) and 
$, we are effectively assuming that $ evolves very slowly 
over the time interval between snapshots. This is a reason- 
able assumption on large scales, since, as discussed in i|2]the 
time derivative is close to zero for most of the evolution of 
the Universe and only weakly evolving away from this at 
later times. On smaller scales this may be a less reasonable 
approximation, however, we are still using the fully nonlin- 
ear gravitational potential field. 
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